library(tidyverse)
library(phyloseq)
library(speedyseq)
library(ggrepel)
library(ampvis2)
library(plotly)
library(microbiome)
options(getClass.msg=FALSE) # https://github.com/epurdom/clusterExperiment/issues/66
#this fixes an error message that pops up because the class 'Annotated' is defined in two different packages
'%!in%' <- function(x,y)!('%in%'(x,y))
source("https://raw.githubusercontent.com/fconstancias/DivComAnalyses/master/R/phyloseq_taxa_tests.R")
source("https://raw.githubusercontent.com/fconstancias/DivComAnalyses/master/R/phyloseq_normalisation.R")
## Loading required package: scales
##
## Attaching package: 'scales'
## The following object is masked from 'package:microbiome':
##
## alpha
## The following object is masked from 'package:purrr':
##
## discard
## The following object is masked from 'package:readr':
##
## col_factor
## Loading required package: reshape2
##
## Attaching package: 'reshape2'
## The following object is masked from 'package:tidyr':
##
## smiths
source("https://raw.githubusercontent.com/fconstancias/DivComAnalyses/master/R/phyloseq_alpha.R")
source("https://raw.githubusercontent.com/fconstancias/DivComAnalyses/master/R/phyloseq_beta.R")
source("https://raw.githubusercontent.com/fconstancias/DivComAnalyses/master/R/phyloseq_heatmap.R")
ps = "data/processed/physeq_update_11_1_21.RDS"
ps %>%
here::here() %>%
readRDS() %>%
phyloseq_get_strains_fast() %>%
phyloseq_remove_chloro_mitho() -> physeq
## Joining, by = "ASV"
physeq %>%
subset_samples(Experiment == "Continuous") %>%
subset_samples(Paul %!in% c("Paul")) %>%
subset_samples(Reactor != "IR2") -> ps_PolyFermS
We will be analysing only the PolyFermS samples here so take a subset of the physeq object.
physeq %>%
subset_samples(Experiment == "Continuous") %>%
subset_samples(Paul %!in% c("Paul")) %>%
subset_samples(Reactor != "IR2") -> ps_polyFermS
sample_data(ps_polyFermS)$Reactor <- fct_relevel(sample_data(ps_polyFermS)$Reactor, "IR1", "CR", "TR1", "TR2","TR3", "TR4", "TR5", "TR6")
sample_data(ps_polyFermS)$Treatment <- fct_relevel(sample_data(ps_polyFermS)$Treatment, "UNTREATED", "CTX+HV292.1", "CTX","HV292.1","VAN+CCUG59168", "VAN", "CCUG59168")
sample_data(ps_polyFermS)$Reactor_Treatment <- fct_relevel(sample_data(ps_polyFermS)$Reactor_Treatment, "IR1_UNTREATED","CR_UNTREATED", "CR_CTX", "CR_VAN", "TR1_CTX+HV292.1","TR2_CTX", "TR3_HV292.1", "TR5_VAN+CCUG59168", "TR4_VAN", "TR6_CCUG59168")
ps_polyFermS %>%
rarefy_even_depth(sample.size = 4576,
rngseed = 123) -> ps_rare
## `set.seed(123)` was used to initialize repeatable random subsampling.
## Please record this for your records so others can reproduce.
## Try `set.seed(123); .Random.seed` for the full vector
## ...
## 16 samples removedbecause they contained fewer reads than `sample.size`.
## Up to first five removed samples are:
## CR-40-S166ETR1-30-S178ETR1-42-S194ETR2-30-S195IR1-40-S197
## ...
## 50OTUs were removed because they are no longer
## present in any sample after random subsampling
## ...
ps_rare %>%
subset_samples(Enrichment == "NotEnriched") %>%
subset_samples(Sample_description %!in% c("TR5-15", "TR4-1")) %>% # heatmap reavels they might be comming from VAN treated or enrichments
subset_samples(Day_from_Inoculum >= 30 & Day_from_Inoculum <= 38) -> ps_polyFermS_rare
"data/raw/hplc Fermentation (Salvato automaticamente).xlsx" %>%
readxl::read_xlsx(sheet = "All total") -> metabolites
ps_polyFermS_rare@sam_data %>%
data.frame() %>%
rownames_to_column('id') %>%
left_join(
metabolites,
by = c("Sample_description" = "Sample_Id"),
suffix = c(".x", "")) %>%
column_to_rownames('id') %>%
sample_data() -> ps_polyFermS_rare@sam_data
ps_polyFermS_rare
## phyloseq-class experiment-level object
## otu_table() OTU Table: [ 296 taxa and 51 samples ]:
## sample_data() Sample Data: [ 51 samples by 63 sample variables ]:
## tax_table() Taxonomy Table: [ 296 taxa by 8 taxonomic ranks ]:
## phy_tree() Phylogenetic Tree: [ 296 tips and 294 internal nodes ]:
## refseq() DNAStringSet: [ 296 reference sequences ]
## taxa are rows
ps_polyFermS_rare %>%
sample_data() %>%
data.frame() -> df
measures = df %>% dplyr::select(ends_with("mM")) %>% colnames()
df %>%
plot_alphas(measure = measures,
x_group = "Reactor_Treatment",
colour_group = "Reactor",
fill_group = "Reactor",
shape_group = NULL,
facet_group = "Reactor_Treatment",
test_group = "Reactor_Treatment",
test_group_2 = NULL) -> out
out$plot +
facet_null() +
facet_grid(alphadiversiy ~ ., scales = "free") +
ggpubr::rotate_x_text(60) +
scale_fill_viridis_d() +
scale_color_viridis_d() +
ylab("SCFA concentration [mM]")
out$stat %>%
# dplyr::filter(signif == "SIGN") %>%
DT::datatable()
plot_time <- function(df,
measure,
x = "Day_from_Inoculum",
y = "value",
shape = "neg",
fill = "Reactor_Treatment",
group = "Reactor_Treatment",
facet)
{
df %>%
dplyr::filter(alphadiversiy %in% measure) %>%
dplyr::mutate(alphadiversiy = fct_reorder(alphadiversiy, value, .desc = TRUE)) %>%
dplyr::mutate(neg = ifelse(value == 0, "neg", "pos")) %>%
arrange(Day_from_Inoculum) %>%
ggplot(aes_string(x = x,
y = y)) +
geom_jitter(size=0.5, alpha=0.9, aes_string(color = fill, fill = fill, shape = shape), show.legend = TRUE) +
geom_path(inherit.aes = TRUE, aes_string(fill = fill, color = fill, show.legend = FALSE),
size = 0.005,
linetype = "dashed") +
facet_grid(as.formula(facet), scales = "free") +
# geom_vline(xintercept = c(23,39),
# color="black", alpha=0.4) +
geom_smooth(show.legend = FALSE, level = 0.95, alpha=0.05, size = 0.5 ,aes_string(color = fill, fill = fill)) +
# scale_x_continuous(breaks=seq(0,90,10)) +
# scale_y_continuous(labels = scientific,
# limits=c(1e+10, 1e+11), breaks = seq(1e+10, 1e+11, by = 1e+10),
# trans = "log10") +
theme_light() +
scale_color_viridis_d(na.value = "black") +
scale_fill_viridis_d(na.value = "black") -> plot
return(plot + theme(legend.position = "bottom"))
}
out$plot$data %>%
plot_time(measure = c("Total_SCFA_mM", "Acetat_mM", "Butyrat_mM", "Propionat_mM", "Isobutyrat_mM", "Valerat_mM", "Isovalerat_mM", "Succinat_mM"),
facet = c("alphadiversiy ~ ."), shape = NULL) +
labs(x="Day (from Inoculum)", y= "SCFA concentration [mM]",
col=NULL, fill = NULL, shape = NULL) +
scale_shape_manual(values=c(4, 19)) +
geom_vline(xintercept = c(39)) -> p4
p4
p4 %>% ggplotly()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
p4 + geom_boxplot(aes(group = Reactor_Treatment,
color = Reactor_Treatment,
fill = Reactor_Treatment),
alpha = 0.2) -> p5
p5
p4 +
facet_null() +
facet_grid(alphadiversiy ~ Reactor_Treatment, scales = "free") +
scale_color_manual(values = rep("black",8)) +
scale_fill_manual(values = rep("black",8))
df %>%
dplyr::select("Acetat_mM", "Butyrat_mM", "Propionat_mM", "Isobutyrat_mM", "Valerat_mM", "Isovalerat_mM", "Succinat_mM") %>%
# dplyr::select(ends_with("mM") | "Total_SCFA_mM") %>%
drop_na() %>%
# t() %>%
scale(center = TRUE,
scale = TRUE) %>%
dist(method= "euc") -> euc_met
plot_ordination(ps_polyFermS,
ordination = phyloseq::ordinate(ps_polyFermS,
distance = euc_met,
method = "PCoA")) -> pca
pca$layers[[1]] = NULL
pca +
geom_point(size=2,
aes(color = Reactor_Treatment,
fill = NULL,
shape = NULL,
alpha = Day_from_Inoculum)) +
theme_light() +
geom_path(arrow = arrow(type = "open", angle = 30, length = unit(0.15, "inches")),
size = 0.08, linetype = "dashed", inherit.aes = TRUE, aes(group=Reactor_Treatment, color = Reactor_Treatment)) +
scale_alpha_continuous(range=c( 0.9, 0.3)) +
scale_color_viridis_d(na.value = "red") +
scale_fill_viridis_d(na.value = "red") +
scale_shape_manual(values = c(8, 21, 22, 23, 24, 16, 15, 18, 17)) +
theme_classic() -> p5
p5
# Alpha div:
Compute alpha-div metrics:
ps_polyFermS_rare %>%
phyloseq_alphas(phylo = TRUE) -> alpha_df
measures = c("Observed", "diversity_shannon", "evenness_pielou","SES.MPD")
alpha_df %>%
plot_alphas(measure = measures,
x_group = "Reactor_Treatment",
colour_group = "Reactor",
fill_group = "Reactor",
shape_group = NULL,
facet_group = "Reactor_Treatment",
test_group = "Reactor_Treatment",
test_group_2 = NULL) -> out
out$plot +
facet_null() +
facet_grid(alphadiversiy ~ ., scales = "free") +
ggpubr::rotate_x_text(60) +
scale_fill_viridis_d() +
scale_color_viridis_d()
out$stat %>%
# dplyr::filter(signif == "SIGN") %>%
DT::datatable()
out$plot$data %>%
plot_time(measure = measures,
facet = c("alphadiversiy ~ ."), shape = NULL) +
labs(x="Day (from Inoculum)", y= "Alpha-diversity",
col=NULL, fill = NULL, shape = NULL) +
scale_shape_manual(values=c(4, 19)) +
geom_vline(xintercept = c(39)) -> p4
p4
p4 %>% ggplotly()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : span too small. fewer data values than degrees of freedom.
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : pseudoinverse used at 30.965
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : neighborhood radius 2.035
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : reciprocal condition number 0
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : There are other near singularities as well. 25.351
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : span too small. fewer
## data values than degrees of freedom.
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : pseudoinverse used at
## 30.965
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : neighborhood radius
## 2.035
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : reciprocal condition
## number 0
## Warning in predLoess(object$y, object$x, newx = if
## (is.null(newdata)) object$x else if (is.data.frame(newdata))
## as.matrix(model.frame(delete.response(terms(object)), : There are other near
## singularities as well. 25.351
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in simpleLoess(y, x, w, span, degree = degree, parametric =
## parametric, : Chernobyl! trL>n 6
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning in stats::qt(level/2 + 0.5, pred$df): NaNs produced
p4 + geom_boxplot(aes(group = Reactor_Treatment,
color = Reactor_Treatment,
fill = Reactor_Treatment),
alpha = 0.2) -> p5
p5
Compute beta-div metrics:
ps_polyFermS_rare %>%
phyloseq_plot_bdiv(bdiv_list,
m = "CoDa",
seed = 123) -> coda
coda$PCA$layers[[1]] = NULL
coda$PCA + geom_point(size=2,
aes(color = Reactor_Treatment,
fill = Reactor_Treatment,
shape = NULL,
alpha = Day_from_Inoculum)) +
theme_light() +
# geom_path(arrow = arrow(type = "open", angle = 30, length = unit(0.15, "inches")),
# size = 0.05, linetype = "dashed", inherit.aes = TRUE, aes(group=Reactor_Treatment, color = Reactor_Treatment), show.legend = FALSE) +
scale_alpha_continuous(range=c( 0.9, 0.3)) + scale_color_viridis_d(na.value = "black") +
scale_fill_viridis_d(na.value = "black") +
# scale_shape_manual(values = c(8, 21, 22, 23, 24, 16, 15, 18, 17)) +
theme_classic() +
labs(col=NULL, fill = NULL, shape = NULL) + guides(shape=FALSE) -> p1
p1
p1 %>%
ggplotly() -> p1ply
p1ply
ps_polyFermS_rare %>%
phyloseq_compute_bdiv(norm = "pc",
phylo = TRUE,
seed = 123) -> bdiv_list
## Loading required package: ape
## Loading required package: GUniFrac
## Loading required package: vegan
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.5-7
##
## Attaching package: 'vegan'
## The following object is masked from 'package:microbiome':
##
## diversity
## Loading required package: matrixStats
##
## Attaching package: 'matrixStats'
## The following object is masked from 'package:dplyr':
##
## count
bdiv_list$coda <- coda$physeq_clr %>% phyloseq::distance(method = "euclidean")
phyloseq_plot_bdiv(dlist = bdiv_list, # list of distance computed from a phyloseq object
ps_rare = ps_polyFermS_rare, # phyloseq object
m = "PCoA", # PCoA or NMDS
seed = 123, # for reproducibility
axis1 = 1, # axis to plot
axis2 = 2) -> plot_list
## [1] "bray"
## [1] "bjaccard"
## [1] "wjaccard"
## [1] "uunifrac"
## [1] "wunifrac"
## [1] "d_0"
## [1] "d_0.5"
## [1] "coda"
plot_list$bray = NULL; plot_list$d_0 = NULL; plot_list$d_0.5 = NULL
plot_list %>%
phyloseq_plot_ordinations_facet(color_group = "Reactor_Treatment",
shape_group = NULL) +
scale_color_viridis_d(na.value = "black") +
scale_fill_viridis_d(na.value = "black")
plot_list %>%
phyloseq_ordinations_expl_var()
## New names:
## * .id -> .id...1
## * V1 -> V1...2
## * .id -> .id...3
## * V1 -> V1...4
## .id...1 V1...2 .id...3 V1...4
## 1 bjaccard Axis.1 [14.2%] bjaccard Axis.2 [11%]
## 2 wjaccard Axis.1 [23.1%] wjaccard Axis.2 [20.2%]
## 3 uunifrac Axis.1 [15.3%] uunifrac Axis.2 [12.4%]
## 4 wunifrac Axis.1 [71.6%] wunifrac Axis.2 [12.3%]
## 5 coda Axis.1 [16.2%] coda Axis.2 [12.6%]
ps_polyFermS_rare %>%
phyloseq_distance_boxplot(bdiv_list$coda ,
d = "Reactor_Treatment") -> dist_box
dist_box$plot
lapply(
bdiv_list,
FUN = physeq_pairwise_permanovas,
physeq = ps_polyFermS_rare,
compare_header = "Reactor_Treatment",
n_perm = 999,
strat = FALSE
) %>%
bind_rows(.id = "Distance") %>%
mutate_if(is.numeric, round, 3) %>%
# filter(! terms %in% (c("Residuals", "Total"))) %>%
DT::datatable()
ps_polyFermS_rare %>%
# subset_samples(Experiment == "CCUG59168" | Experiment == "HV292.1" | Experiment == "Cecum" | Experiment == "Continuous" & Reactor %in% c("IR1", "CR")) %>%
# rarefy_even_depth(sample.size = 2574,rngseed = 123) %>%
phyloseq_ampvis_heatmap(transform = "percent",
group_by = "Day_from_Inoculum",
facet_by = c("Enrichment", "Phase", "Reactor", "Treatment", "Experiment", "Reactor_Treatment" ),
tax_aggregate = "Class",
tax_add = NULL,
ntax = 60) -> p
## Warning: There are only 7 taxa, showing all
p + facet_grid( ~ Reactor_Treatment + Phase , scales = "free", space = "free") +
scale_fill_viridis_c(breaks = c(0, 0.01, 1, 10, 50, 75, 100),
labels = c(0, 0.01, 1, 10, 50, 75, 100),
trans = scales::pseudo_log_trans(sigma = 0.001),
na.value = 'transparent') -> p2
## Scale for 'fill' is already present. Adding another scale for 'fill', which
## will replace the existing scale.
p2
ps_polyFermS_rare %>%
# subset_samples(Experiment == "CCUG59168" | Experiment == "HV292.1" | Experiment == "Cecum" | Experiment == "Continuous" & Reactor %in% c("IR1", "CR")) %>%
# rarefy_even_depth(sample.size = 2574,rngseed = 123) %>%
phyloseq_ampvis_heatmap(transform = "percent",
group_by = "Day_from_Inoculum",
facet_by = c("Enrichment", "Phase", "Reactor", "Treatment", "Experiment", "Reactor_Treatment" ),
tax_aggregate = "Genus",
tax_add = NULL,
ntax = 60) -> p
p + facet_grid( ~ Reactor_Treatment + Phase , scales = "free", space = "free") +
scale_fill_viridis_c(breaks = c(0, 0.01, 1, 10, 50, 75, 100),
labels = c(0, 0.01, 1, 10, 50, 75, 100),
trans = scales::pseudo_log_trans(sigma = 0.001),
na.value = 'transparent') -> p2
## Scale for 'fill' is already present. Adding another scale for 'fill', which
## will replace the existing scale.
p2
ps_polyFermS_rare %>%
# subset_samples(Experiment == "CCUG59168" | Experiment == "HV292.1" | Experiment == "Cecum" | Experiment == "Continuous" & Reactor %in% c("IR1", "CR")) %>%
# rarefy_even_depth(sample.size = 2574,rngseed = 123) %>%
phyloseq_ampvis_heatmap(transform = "percent",
group_by = "Day_from_Inoculum",
facet_by = c("Enrichment", "Phase", "Reactor", "Treatment", "Experiment", "Reactor_Treatment" ),
tax_aggregate = "Species",
tax_add = NULL,
ntax = 60) -> p
p + facet_grid( ~ Reactor_Treatment + Phase , scales = "free", space = "free") +
scale_fill_viridis_c(breaks = c(0, 0.01, 1, 10, 50, 75, 100),
labels = c(0, 0.01, 1, 10, 50, 75, 100),
trans = scales::pseudo_log_trans(sigma = 0.001),
na.value = 'transparent') -> p2
## Scale for 'fill' is already present. Adding another scale for 'fill', which
## will replace the existing scale.
p2
paste0(here::here(),
"/data/processed/",
"stab",
"_",
format(Sys.time(), "%Y%b%d")
,".RData") %>% save.image()
# load("/Users/fconstan/Projects/EZe/ASV/260420.RData")
sessionInfo()
## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Mojave 10.14.6
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] GUniFrac_1.1 matrixStats_0.57.0 vegan_2.5-7
## [4] lattice_0.20-41 permute_0.9-5 ape_5.4-1
## [7] metagMisc_0.0.4 reshape2_1.4.4 scales_1.1.1
## [10] microbiome_1.10.0 plotly_4.9.3 ampvis2_2.6.4
## [13] ggrepel_0.9.0 speedyseq_0.5.3.9001 phyloseq_1.32.0
## [16] forcats_0.5.0 stringr_1.4.0 dplyr_1.0.3
## [19] purrr_0.3.4 readr_1.4.0 tidyr_1.1.2
## [22] tibble_3.0.5 ggplot2_3.3.3 tidyverse_1.3.0.9000
##
## loaded via a namespace (and not attached):
## [1] readxl_1.3.1 backports_1.2.1 plyr_1.8.6
## [4] igraph_1.2.6 lazyeval_0.2.2 splines_4.0.3
## [7] crosstalk_1.1.1 digest_0.6.27 foreach_1.5.1
## [10] htmltools_0.5.1 fansi_0.4.2 magrittr_2.0.1
## [13] cluster_2.1.0 doParallel_1.0.16 openxlsx_4.2.3
## [16] Biostrings_2.56.0 modelr_0.1.8 prettyunits_1.1.1
## [19] colorspace_2.0-0 rvest_0.3.6 haven_2.3.1
## [22] xfun_0.20 crayon_1.3.4 jsonlite_1.7.2
## [25] survival_3.2-7 iterators_1.0.13 glue_1.4.2
## [28] gtable_0.3.0 zlibbioc_1.34.0 XVector_0.28.0
## [31] car_3.0-10 Rhdf5lib_1.10.0 BiocGenerics_0.34.0
## [34] abind_1.4-5 DBI_1.1.0 rstatix_0.6.0
## [37] Rcpp_1.0.6 viridisLite_0.3.0 progress_1.2.2
## [40] foreign_0.8-81 stats4_4.0.3 DT_0.17
## [43] truncnorm_1.0-8 htmlwidgets_1.5.3 httr_1.4.2
## [46] RColorBrewer_1.1-2 ellipsis_0.3.1 pkgconfig_2.0.3
## [49] NADA_1.6-1.1 farver_2.0.3 dbplyr_2.0.0
## [52] here_1.0.1 tidyselect_1.1.0 labeling_0.4.2
## [55] rlang_0.4.10 munsell_0.5.0 cellranger_1.1.0
## [58] tools_4.0.3 cli_2.2.0 generics_0.1.0
## [61] ade4_1.7-16 broom_0.7.3 evaluate_0.14
## [64] biomformat_1.16.0 PhyloMeasures_2.1 yaml_2.2.1
## [67] knitr_1.30 fs_1.5.0 zip_2.1.1
## [70] nlme_3.1-151 xml2_1.3.2 compiler_4.0.3
## [73] rstudioapi_0.13 curl_4.3 ggsignif_0.6.0
## [76] zCompositions_1.3.4 reprex_0.3.0 stringi_1.5.3
## [79] Matrix_1.3-2 multtest_2.44.0 vctrs_0.3.6
## [82] pillar_1.4.7 lifecycle_0.2.0 data.table_1.13.6
## [85] patchwork_1.1.1 R6_2.5.0 network_1.16.1
## [88] rio_0.5.16 IRanges_2.22.2 codetools_0.2-18
## [91] ggnet_0.1.0 MASS_7.3-53 assertthat_0.2.1
## [94] picante_1.8.2 rhdf5_2.32.2 rprojroot_2.0.2
## [97] withr_2.4.0 S4Vectors_0.26.1 mgcv_1.8-33
## [100] parallel_4.0.3 hms_0.5.3 grid_4.0.3
## [103] rmarkdown_2.6 carData_3.0-4 Rtsne_0.15
## [106] ggpubr_0.4.0 Biobase_2.48.0 lubridate_1.7.9.2